home *** CD-ROM | disk | FTP | other *** search
/ CDUTIL 13 / CDUTIL #13 Julio 1995.iso / windows / acadcom / ads / sample / magnets.c < prev    next >
Encoding:
Text File  |  1995-02-08  |  57.1 KB  |  1,610 lines

  1. /* Next available MSG number is 114 */
  2.  
  3. /*    
  4.  
  5.    MAGNETS.C
  6.  
  7.    Copyright (C) 1989, 1990, 1991, 1992, 1993, 1994 by Autodesk, Inc.
  8.  
  9.    Permission to use, copy, modify, and distribute this software in 
  10.    object code form for any purpose and without fee is hereby granted, 
  11.    provided that the above copyright notice appears in all copies and 
  12.    that both that copyright notice and the limited warranty and 
  13.    restricted rights notice below appear in all supporting 
  14.    documentation.
  15.  
  16.    AUTODESK PROVIDES THIS PROGRAM "AS IS" AND WITH ALL FAULTS.  
  17.    AUTODESK SPECIFICALLY DISCLAIMS ANY IMPLIED WARRANTY OF 
  18.    MERCHANTABILITY OR FITNESS FOR A PARTICULAR USE.  AUTODESK, INC.
  19.    DOES NOT WARRANT THAT THE OPERATION OF THE PROGRAM WILL BE 
  20.    UNINTERRUPTED OR ERROR FREE.
  21.  
  22.    Use, duplication, or disclosure by the U.S. Government is subject to 
  23.    restrictions set forth in FAR 52.227-19 (Commercial Computer 
  24.    Software - Restricted Rights) and DFAR 252.227-7013(c)(1)(ii) 
  25.    (Rights in Technical Data and Computer Software), as applicable.
  26.     
  27.    .
  28.  
  29.       DESCRIPTION:
  30.  
  31.  
  32.                   S Q U I R M I N G   M A G N E T S   !!!
  33.  
  34.         Magnetically perturbed pendulum chaos demonstrator for AutoCAD
  35.  
  36.         A sample ADS application.
  37.  
  38.         Designed and implemented by John Walker in October of 1989.
  39.  
  40.         The  physics  of  the  simulation  is  abstract,  not precise.
  41.         First, I use arbitrary units for the mass and strength of  the
  42.         magnets,  chosen so that values in the range from -10 to 10 do
  43.         reasonable things.  Second,  I  treat  the  magnets  as  point
  44.         sources  of  magnetic  field  lines--as  if they were magnetic
  45.         monopoles (or  the  whole  simulation  were  of  electrostatic
  46.         charges  rather than magnets).  Since real magnets always have
  47.         two poles, and since magnets are of substantial size  compared
  48.         to  the  scale of the simulation, there are substantial dipole
  49.         and quadupole moments to the force exerted on the pendulum  in
  50.         the actual apparatus.  None of these considerations are really
  51.         relevant to any of the phenomena  of  attractors,  orbits,  or
  52.         chaotic  processes,  only to the fidelity of the simulation to
  53.         an actual pendulum.
  54.  
  55.         The following commands are implemented in this file:
  56.  
  57.         DEMO        Creates  one  of  a series of standard demo models
  58.                     when  executed  in  a  new,  empty drawing.  These
  59.                     models may be modified with the MAGNET and MAGEDIT
  60.                     commands just like models entered manually.
  61.  
  62.         MAGNET      Creates  a  new  magnet.   You're invited to enter
  63.                     its   position   (by   any  means  of  co-ordinate
  64.                     specification),  and  its  strength  in  arbitrary
  65.                     units.  If the strength is positive, it repels the
  66.                     pendulum; if positive, it attracts it.
  67.  
  68.         MAGEDIT     Lets  you  modify  the  strength  of  an  existing
  69.                     magnet.   Pick  a  single  magnet  by pointing.  A
  70.                     dialogue is displayed with the properties  of  the
  71.                     magnet.   Change them as you wish, then pick OK to
  72.                     update  the  properties  of  the  magnet  in   the
  73.                     database.   If  you pick Cancel, the magnet is not
  74.                     modified.
  75.  
  76.         RESET       When  a  simulation  is  run,  it  halts after the
  77.                     specified number of steps or  simulated  time  has
  78.                     elapsed.    A   subsequent  RUN  command  normally
  79.                     resumes the simulation from the point at which the
  80.                     last  stopped.  RESET erases all motion paths from
  81.                     the screen and sets the  simulation  back  to  the
  82.                     start.    RUN   after   a  RESET  will  begin  the
  83.                     simulation from the initial state.
  84.  
  85.         RUN         Starts the simulation.  You're  asked  to  specify
  86.                     the length of the simulation as either a number of
  87.                     motion   steps  or  by  the  number  of  simulated
  88.                     seconds.  If you enter  a  positive  number,  it's
  89.                     taken  as  a  step count.  Negative numbers (which
  90.                     may include decimal fractions)  specify  simulated
  91.                     time  in  seconds.  You can terminate a simulation
  92.                     with the Control C key.  RUN can also  be  invoked
  93.                     as  a  function,  (C:RUN <length>), where <length>
  94.                     specifies the simulation length by a  positive  or
  95.                     negative  number  as described above.  When called
  96.                     as a function, C:RUN returns  the  simulation  end
  97.                     time as its result.
  98.  
  99.         SETMAG      The  magnets  simulator has several variables that
  100.                     control  its operation.  These variables are saved
  101.                     with the drawing and may be inspected and modified
  102.                     with   the  SETMAG  command.   SETMAG  displays  a
  103.                     dialogue in  which  you  may  change  any  of  the
  104.                     following variables:
  105.  
  106.                        Output to display only?
  107.                           This is set to Yes or  No  (True/False,  and
  108.                           1/0  are  also  accepted).   If  set  to the
  109.                           default value of Yes, motion paths are drawn
  110.                           using  temporary  vectors within the current
  111.                           viewport.  If the picture is  REDRAWn,  they
  112.                           disappear.   If  set to No, motion paths are
  113.                           added to the drawing as LINE entities.  This
  114.                           causes paths to appear in all viewports, and
  115.                           you can  ZOOM  on  sections  of  a  path  to
  116.                           examine   it   in   greater  detail.   Paths
  117.                           represented as  LINEs  are  saved  with  the
  118.                           drawing.  Generating LINE entities for paths
  119.                           is much slower than just drawing them on the
  120.                           screen,  so  choose  this mode only when you
  121.                           need it.
  122.  
  123.                        Display step number?
  124.                           If this mode is  "Yes"  (the  default),  the
  125.                           simulation  step  number is displayed in the
  126.                           coordinate status  line  as  the  simulation
  127.                           progresses.
  128.  
  129.                        Display time?
  130.                           If "Yes" (the default) the simulated time in
  131.                           years is displayed in the coordinate  status
  132.                           line.
  133.  
  134.                        Step size?
  135.                           This  real  number specifies the factor used
  136.                           to determine the  size  of  the  integration
  137.                           steps  used in calculating the motion of the
  138.                           magnets.  The  time  step  is  specified  in
  139.                           terms  of seconds.  The default value of 0.1
  140.                           works  well  for   reasonably   well-behaved
  141.                           simulations.   If you find that a simulation
  142.                           is taking too long, you can speed it  up  by
  143.                           increasing  the step size, but be aware that
  144.                           the results you see may  not  be  physically
  145.                           accurate.     Increasing   the   step   size
  146.                           magnifies the  inaccuracies  of  modeling  a
  147.                           continuous    process   such   as   magnetic
  148.                           attraction by discrete steps.   In  general,
  149.                           you  can  increase  the step size as long as
  150.                           you  obtain  the  same  results.   When  the
  151.                           outcome  begins to vary, you've set the step
  152.                           size too large.
  153.  
  154.                        Dual pendulum offset?
  155.                           To  observe  the  sensitive  dependence   on
  156.                           initial conditions that characterises chaos,
  157.                           you can simulate two pendulums  which  start
  158.                           with  a  small  initial  displacement.   The
  159.                           pendulums do not interact with one  another,
  160.                           but  respond  to  the  magnets placed on the
  161.                           baseplate.  If the offset is zero, a  single
  162.                           pendulum  is  simulated, resulting in faster
  163.                           execution.
  164. */
  165.  
  166. #include   <stdio.h>
  167. #include   <string.h>
  168. #include   <ctype.h>
  169. #include   <math.h>
  170. #include    <assert.h>
  171.  
  172. #include   "adslib.h"
  173.  
  174. /*  Standard drawing object names  */
  175.  
  176. /* Utility frozen layer for information */
  177. #define FrozenLayer /*MSG1*/"FROZEN-SOLID"
  178. #define OrbitLayer  /*MSG2*/"PENDPATH"     /* Layer for pendulum path traces */
  179.  
  180. /*  Data types  */
  181.  
  182. typedef enum {False = 0, True = 1} Boolean;
  183.  
  184. #define HANDLEN  18                   /* String long enough to hold a handle */
  185.  
  186. /* Definitions  to  wrap  around  submission  of  AutoCAD commands to
  187.    prevent their being echoed.  */
  188.  
  189. #define Cmdecho  False                /* Make True for debug command output */
  190.  
  191. #define CommandB()  { struct resbuf rBc, rBb, rBu, rBh; \
  192.         ads_getvar(/*MSG0*/"CMDECHO", &rBc); \
  193.         ads_getvar(/*MSG0*/"BLIPMODE", &rBb); \
  194.         ads_getvar(/*MSG0*/"HIGHLIGHT", &rBh); \
  195.         rBu.restype = RTSHORT; \
  196.         rBu.resval.rint = (int) Cmdecho; \
  197.         ads_setvar(/*MSG0*/"CMDECHO", &rBu); \
  198.         rBu.resval.rint = (int) False; \
  199.         ads_setvar(/*MSG0*/"BLIPMODE", &rBu); \
  200.         ads_setvar(/*MSG0*/"HIGHLIGHT", &rBu)
  201.  
  202. #define CommandE()  ads_setvar(/*MSG0*/"CMDECHO", &rBc); \
  203.                     ads_setvar(/*MSG0*/"BLIPMODE", &rBb); \
  204.                     ads_setvar(/*MSG0*/"HIGHLIGHT", &rBh); }
  205.  
  206. /*  Definitions  that permit you to push and pop system variables with
  207.     minimal complexity.  These don't work  (and  will  cause  horrible
  208.     crashes  if  used)  with  string  variables,  but since all string
  209.     variables are read-only, they cannot be saved and restored in  any
  210.     case.  */
  211.  
  212. #define PushVar(var, cell, newval, newtype) { struct resbuf cell, cNeW; \
  213.         ads_getvar(var, &cell); cNeW.restype = cell.restype;             \
  214.         cNeW.resval.newtype = newval; ads_setvar(var, &cNeW)
  215.  
  216. #define PopVar(var, cell) ads_setvar(var, &cell); }
  217.  
  218. /* Set point variable from three co-ordinates */
  219.  
  220. #define Spoint(pt, x, y, z)  pt[X] = (x);  pt[Y] = (y);  pt[Z] = (z)
  221.  
  222. /* Copy point  */
  223.  
  224. #define Cpoint(d, s)   d[X] = s[X];  d[Y] = s[Y];  d[Z] = s[Z]
  225.  
  226. /* Generation parameters for objects */
  227.  
  228. #define PendLength  20.0              /* Pendulum length */
  229.  
  230. /* Utility definition to get an  array's  element  count  (at  compile
  231.    time).   For  example:  
  232.  
  233.        int  arr[] = {1,2,3,4,5};
  234.        ... 
  235.        printf("%d", ELEMENTS(arr));
  236.  
  237.    would print a five.  ELEMENTS("abc") can also be used to  tell  how
  238.    many  bytes are in a string constant INCLUDING THE TRAILING NULL. */
  239.  
  240. #define ELEMENTS(array) (sizeof(array)/sizeof((array)[0]))
  241.  
  242. /* Utility definitions */
  243.  
  244. #ifndef abs
  245. #define abs(x)      ((x)<0 ? -(x) : (x))
  246. #endif  /* abs */
  247.  
  248. /*  Forward functions  */
  249.  
  250. void    main _((int, char **));
  251. Boolean funcload _((void));
  252. char    *alloc _((unsigned));
  253. void    defmagblk _((void));
  254. struct resbuf *entitem _((struct resbuf *, int));
  255. void    entinfo _((ads_name, char *, ads_point, ads_real *, int *));
  256. void    partext _((void));
  257. void    sumvec _((ads_real *, ads_real *, ads_real, ads_real *));
  258. void    varblockdef _((void));
  259. void    savemodes _((void));
  260. Boolean boolvar _((char *, char *));
  261. void    varset _((void));
  262. Boolean initacad _((Boolean));
  263. void    addmag _((char *, ads_point, ads_real));
  264. void    magnet _((void));
  265. void    demo _((void));
  266. void    magedit _((void));
  267. void    reset _((void));
  268. void    run _((void));
  269. void    setmag _((void));
  270.  
  271. static Boolean functional;            /* C:command is returning result */
  272. static ads_real numsteps = 1000;      /* Default number of steps to run */
  273. static double simtime = 0.0;          /* Simulated time */
  274. static long stepno = 0;               /* Step number */
  275.  
  276. /*  Command definition and dispatch table.  */
  277.  
  278. struct {
  279.         char *cmdname;
  280.         void (*cmdfunc)();
  281. } cmdtab[] = {
  282. /*        Name         Function  */
  283. {/*MSG3*/"DEMO",       demo},         /* Create demo case */
  284. {/*MSG4*/"MAGNET",     magnet},       /* Create new magnet */
  285. {/*MSG5*/"MAGEDIT",    magedit},      /* Edit existing magnet */
  286. {/*MSG6*/"RESET",      reset},        /* Erase motion paths */
  287. {/*MSG7*/"RUN",        run},          /* Run simulation */
  288. {/*MSG8*/"SETMAG",     setmag}        /* Set mode variables */
  289. };
  290.  
  291. /*  Particle structure.  */
  292.  
  293. typedef struct {
  294.         char partname[32];            /* Particle name */
  295.         ads_point position;           /* Location in space */
  296.         ads_point velocity;           /* Velocity vector */
  297.         ads_point acceleration;       /* Acceleration vector */
  298.         ads_point lastpos;            /* Last plotted position */
  299.         ads_real strength;            /* Magnetic strength */
  300.         ads_real radius;              /* Radius */
  301.         int colour;                   /* Entity colour */
  302.         char partblock[HANDLEN];      /* Block defining particle */
  303. } particle;
  304.  
  305. static particle pt;                   /* Static particle structure */
  306. static particle *ptab = NULL;         /* Particle table */
  307. static int nptab = 0;                 /* Number of in-memory particles */
  308.  
  309. static int ready = 0;                 /* don't run till ready */
  310.  
  311. /*  Particle definition block attributes.  */
  312.  
  313. #define ParticleBlock  /*MSG9*/"PARTICLE"  /* Particle block name */
  314. struct {                              /* Attribute tag table */
  315.         char *attname;                /* Attribute tag name */
  316.         ads_real *attvar;             /* Variable address */
  317.         char *attprompt;              /* Prompt for attribute */
  318. } partatt[] = {
  319.         {/*MSG10*/"MAGSTRENGTH", &pt.strength, /*MSG11*/"Magnetic strength"}
  320.        };
  321.  
  322. /*  Modal  variables.   Default   initial   values   below   are   for
  323.     documentation  only.   The actual defaults are reset in initacad()
  324.     upon entry to the drawing editor, so that they will be placed with
  325.     the   default  values  in  the  modal  variable  block  if  it  is
  326.     subsequently created for a new drawing.  */
  327.  
  328. static Boolean drawonly = True,       /* DRAWONLY: Draw, but don't add entities
  329.                                                    if 1.  Make LINEs if 0. */
  330.                showstep = True,       /* SHOWSTEP: Show step in status line. */
  331.                showtime = True;       /* SHOWTIME: Show time in status line. */
  332. static ads_real stepsize = 0.1,       /* STEPSIZE: Motion step size factor. */
  333.                 dualpend = 0.0;       /* DUALPEND: Offset of dual pendulum */
  334. static int ntestp = 1;                /* Number of test particles */
  335.  
  336. /*  Modal attribute definition.  */
  337.  
  338. #define ModalBlock  /*MSG12*/"GRAVITY_MODES" /* Mode variable block name */
  339. struct {
  340.         char *attname;                /* Attribute tag name */
  341.         int attvari;                  /* Variable index */
  342.         char *attprompt;              /* Prompt for variable */
  343. } varatt[] = {
  344.         {/*MSG13*/"DRAWONLY", 1, /*MSG14*/"Output to display only? "},
  345.         {/*MSG15*/"SHOWSTEP", 3, /*MSG16*/"Display step number?    "},
  346.         {/*MSG17*/"SHOWTIME", 2, /*MSG18*/"Display time?           "},
  347.         {/*MSG19*/"STEPSIZE", 4, /*MSG20*/"Step size?              "},
  348.         {/*MSG21*/"DUALPEND", 5, /*MSG22*/"Dual pendulum offset?   "}
  349.        };
  350.  
  351. #ifdef ACRXAPP
  352. int acrxEntryPoint(int stat, void* ptr)
  353. {
  354.   int scode;
  355.   int cindex;
  356.   int functional;
  357.         scode = RSRSLT;               /* Default return code */
  358.  
  359.         switch (stat) {
  360.  
  361.         case 3: /*RQXLOAD */           /* Load functions.  Called at the start
  362.                                          of the drawing editor.  Re-initialise
  363.                                          the application here. */
  364.             scode = -(funcload() ? RSRSLT : RSERR);
  365.             break;
  366.  
  367.         case 5: /*RQSUBR:*/             /* Evaluate external lisp function */
  368.             cindex = ads_getfuncode();
  369.             functional = False;
  370.             if (!initacad(False)) {
  371.                 ads_printf(/*MSG24*/"\nUnable to initialise application.\n");
  372.             } else {
  373.  
  374.                 /* Execute the command from the command table with
  375.                    the index associated with this function. */
  376.  
  377.                 if (cindex > 0) {
  378.                     cindex--;
  379.                     assert(cindex < ELEMENTS(cmdtab));
  380.                     (*cmdtab[cindex].cmdfunc)();
  381.                 }
  382.             }
  383.             if (!functional)
  384.                 ads_retvoid();        /* Void result */
  385.             break;
  386.          }
  387.  
  388.   return 0;
  389.   }
  390. #else
  391. /* MAIN -- the main routine */
  392.  
  393. void main(argc, argv)
  394.   int argc;
  395.   char *argv[];
  396. {
  397.     int stat, cindex, scode = RSRSLT;
  398.     char errmsg[80];
  399.  
  400.     ads_init(argc, argv);             /* Initialise the application */
  401.  
  402.     /* Main dispatch loop. */
  403.  
  404.     while (True) {
  405.  
  406.         if ((stat = ads_link(scode)) < 0) {
  407.             sprintf(errmsg,
  408.                     /*MSG23*/"MAGNETS: bad status from ads_link() = %d\n",
  409.                     stat);
  410. #ifdef Macintosh
  411.             macalert(errmsg);
  412. #else
  413.             puts(errmsg);
  414.             fflush(stdout);
  415. #endif /* Macintosh */
  416.             exit(1);
  417.         }
  418.  
  419.         scode = RSRSLT;               /* Default return code */
  420.  
  421.         switch (stat) {
  422.  
  423.         case RQXLOAD:                 /* Load functions.  Called at the start
  424.                                          of the drawing editor.  Re-initialise
  425.                                          the application here. */
  426.             scode = -(funcload() ? RSRSLT : RSERR);
  427.             break;
  428.  
  429.         case RQSUBR:                  /* Evaluate external lisp function */
  430.             cindex = ads_getfuncode();
  431.             functional = False;
  432.             if (!initacad(False)) {
  433.                 ads_printf(/*MSG24*/"\nUnable to initialise application.\n");
  434.             } else {
  435.  
  436.                 /* Execute the command from the command table with
  437.                    the index associated with this function. */
  438.  
  439.                 if (cindex > 0) {
  440.                     cindex--;
  441.                     assert(cindex < ELEMENTS(cmdtab));
  442.                     (*cmdtab[cindex].cmdfunc)();
  443.                 }
  444.             }
  445.             if (!functional)
  446.                 ads_retvoid();        /* Void result */
  447.             break;
  448.  
  449.         default:
  450.             break;
  451.         }
  452.     }
  453. }
  454. #endif
  455.  
  456. /* FUNCLOAD  --  Load external functions into AutoLISP */
  457.  
  458. static Boolean funcload()
  459. {
  460.     char ccbuf[40];
  461.     int i;
  462.  
  463.     strcpy(ccbuf, /*MSG0*/"C:");
  464.     for (i = 0; i < ELEMENTS(cmdtab); i++) {
  465.         strcpy(ccbuf + 2, cmdtab[i].cmdname);
  466.         ads_defun(ccbuf, i + 1);
  467.     }
  468.  
  469.     return initacad(True);            /* Reset AutoCAD initialisation */
  470. }
  471.  
  472. /*  ALLOC  --  Allocate storage and fail if it can't be obtained.  */
  473.  
  474. static char *alloc(nbytes)
  475.   unsigned nbytes;
  476. {
  477.     char *cp;
  478.  
  479.     if ((cp = malloc(nbytes)) == NULL) {
  480.         ads_printf(/*MSG25*/"Out of memory");
  481.     }
  482.     return cp;
  483. }
  484.  
  485. /*  DEFMAGBLK  --  Create the magnet definition block. */
  486.  
  487. static void defmagblk()
  488. {
  489.     ads_point centre, ucentre;
  490.     ads_real radius = 0.25;
  491.     ads_point ax, ax1;
  492.     struct resbuf rb1, rb2, ocolour;
  493.     ads_name e2;
  494.     int i;
  495.     ads_name matbss;
  496.     ads_name ename;
  497.     ads_point atx;
  498.  
  499.     Spoint(ucentre, 4, 4, 0);
  500.     rb1.restype = rb2.restype = RTSHORT;
  501.     rb1.resval.rint = 1;              /* From UCS */
  502.     rb2.resval.rint = 0;              /* To world */
  503.     ads_trans(ucentre, &rb1, &rb2, False, centre);
  504.     CommandB();
  505.     ads_getvar(/*MSG0*/"CECOLOR", &ocolour);
  506.  
  507.     /* AutoCAD  currently returns  "human readable" colour strings
  508.        like "1 (red)" for the standard colours.  Trim  the  string
  509.        at  the  first space to guarantee we have a valid string to
  510.        restore the colour later.  */
  511.     if (strchr(ocolour.resval.rstring, ' ') != NULL)
  512.         *strchr(ocolour.resval.rstring, ' ') = EOS;
  513.  
  514.     ads_command(RTSTR, /*MSG0*/"_.ELEV", RTREAL, -0.5, RTREAL, 0.5, RTNONE);
  515.     Spoint(ax, 0.0, 0.0, -0.5);
  516.     ads_command(RTSTR, /*MSG0*/"_.COLOUR", RTSTR, /*MSG0*/"_BLUE", RTNONE);
  517. #define PcX   3
  518.     ads_command(RTSTR, /*MSG0*/"_.Circle",
  519.                 RT3DPOINT, ax, RTREAL, PendLength / PcX,
  520.                 RTNONE);
  521.     ads_command(RTSTR, /*MSG0*/"_.ELEV", RTREAL, 0.0, RTREAL, 0.0, RTNONE);
  522. #define PlX   40.0
  523.     Spoint(ax, -PendLength / PlX, -PendLength / PlX, 0.0);
  524.     Spoint(ax1, PendLength / PlX, PendLength / PlX, 0.0);
  525.     ads_command(RTSTR, /*MSG*/"_.COLOUR", RTSTR, /*MSG0*/"_WHITE", RTNONE);
  526.     ads_command(RTSTR, /*MSG*/"_.LINE", RT3DPOINT, ax, RT3DPOINT, ax1,
  527.                 RTSTR, "", RTNONE);
  528.     Spoint(ax, -PendLength / PlX, PendLength / PlX, 0.0);
  529.     Spoint(ax1, PendLength / PlX, -PendLength / PlX, 0.0);
  530.     ads_command(RTSTR, /*MSG0*/"_.LINE", RT3DPOINT, ax, RT3DPOINT, ax1,
  531.                 RTSTR, "", RTNONE);
  532.  
  533.     ads_command(RTSTR, /*MSG0*/"_.COLOUR", RTSTR, /*MSG0*/"_BYBLOCK", RTNONE);
  534.     rb1.resval.rint = 0;              /* From world */
  535.     rb2.resval.rint = 1;              /* To new UCS */
  536.     ads_trans(centre, &rb1, &rb2, False, centre);
  537.     ax[X] = ax1[X] = centre[X];
  538.     ax[Y] = centre[Y] + radius;
  539.     ax[Z] = ax1[Z] = centre[Z];
  540.     ax1[Y] = centre[Y] - radius;
  541.     ads_command(RTSTR, /*MSG0*/"_.ELEV", RTREAL, 0.0, RTREAL, 0.1,
  542.                 RTNONE);
  543.     ads_command(RTSTR, /*MSG0*/"_.Circle", RT3DPOINT, ax, RTREAL, radius,
  544.                 RTNONE);
  545.     ads_command(RTSTR, /*MSG0*/"_.ELEV", RTREAL, 0.0, RTREAL, 0.0, RTNONE);
  546.     ads_command(RTSTR, /*MSG0*/"_.COLOUR",
  547.                 RTSTR, ocolour.resval.rstring, RTNONE);
  548.     free(ocolour.resval.rstring);
  549.     ads_entlast(e2);
  550.  
  551.     /* Create attributes */
  552.  
  553.     Spoint(atx, 4, 2.75, 0);
  554.     ads_ssadd(e2, NULL, matbss);
  555.  
  556.     PushVar(/*MSG0*/"AFLAGS", saflags, 1, rint);  /* Invisible */
  557.     ads_command(RTSTR, /*MSG0*/"_.ATTDEF", RTSTR, "",
  558.                 RTSTR, /*MSG28*/"PARTNAME",
  559.                 RTSTR, /*MSG29*/"Particle name", RTSTR, "", RT3DPOINT, atx,
  560.                 RTREAL, 0.2, RTREAL, 0.0, RTNONE);
  561.     ads_entlast(ename);
  562.     ads_ssadd(ename, matbss, matbss);
  563.  
  564.     for (i = 0; i < ELEMENTS(partatt); i++) {
  565.         atx[Y] -= 0.25;
  566.         ads_command(RTSTR, /*MSG0*/"_.ATTDEF",
  567.                     RTSTR, "", RTSTR, partatt[i].attname,
  568.                     RTSTR, partatt[i].attprompt, RTSTR, "0", RT3DPOINT, atx,
  569.                     RTREAL, 0.2, RTREAL, 0.0, RTNONE);
  570.         ads_entlast(ename);
  571.         ads_ssadd(ename, matbss, matbss);
  572.     }
  573.     PopVar(/*MSG0*/"AFLAGS", saflags);
  574.  
  575.     /* Collect magnet and attributes into a block. */
  576.  
  577.     ads_command(RTSTR, /*MSG0*/"_.BLOCK",
  578.                 RTSTR, ParticleBlock, RT3DPOINT, centre,
  579.                 RTPICKS, matbss, RTSTR, "", RTNONE);
  580.     CommandE();
  581.     ads_ssfree(matbss);
  582. }
  583.  
  584. /*  ENTITEM  --  Search an entity buffer chain and return an item
  585.                  with the specified group code.  */
  586.  
  587. static struct resbuf *entitem(rchain, gcode)
  588.   struct resbuf *rchain;
  589.   int gcode;
  590. {
  591.     while ((rchain != NULL) && (rchain->restype != gcode))
  592.         rchain = rchain->rbnext;
  593.  
  594.     return rchain;
  595. }
  596.  
  597. /*  ENTINFO  --  Obtain information about a particle entity:
  598.  
  599.                        Handle
  600.                        Position
  601.                        Size (from scale factor of unit block)
  602.                        Colour
  603. */
  604.  
  605. static void entinfo(ename, h, p, r, c)
  606.   ads_name ename;
  607.   char *h;
  608.   ads_point p;
  609.   ads_real *r;
  610.   int *c;
  611. {
  612.     struct resbuf *rent, *rh;
  613.  
  614.     rent = ads_entget(ename);
  615.     if ((rh = entitem(rent, 5)) != NULL)
  616.         strcpy(h, rh->resval.rstring);
  617.     else
  618.         h[0] = EOS;
  619.     rh = entitem(rent, 10);
  620.     assert(rh != NULL);
  621.     Cpoint(p, rh->resval.rpoint);
  622.     rh = entitem(rent, 41);
  623.     assert(rh != NULL);
  624.     *r = rh->resval.rreal;
  625.     if ((rh = entitem(rent, 62)) != NULL) {
  626.         *c = rh->resval.rint;
  627.         if (*c == 0)                  /* Naked colour by block? */
  628.             *c = 7;                   /* Forbidden: make it white.  Q.C.D. */
  629.     } else {
  630.         /* This entity derives its colour from the layer.  Get
  631.            the colour from the layer table. */
  632.         *c = 7;
  633.         if ((rh = entitem(rent, 8)) != NULL) {
  634.             if ((rh = ads_tblsearch(/*MSG0*/"LAYER", rh->resval.rstring,
  635.                                     False)) != NULL) {
  636.                 struct resbuf *lh = entitem(rh, 62);
  637.                 if (lh != NULL) {
  638.                     int lc = abs(lh->resval.rint);
  639.                     if (lc > 0 && lc < 256) {
  640.                         *c = lc;
  641.                     }
  642.                 }
  643.                 ads_relrb(rh);
  644.             }
  645.         }
  646.     }
  647.     ads_relrb(rent);
  648. }
  649.  
  650. /*  PARTEXT  --  Extract particles present in drawing and build the
  651.                  in-memory particle table.  */
  652.  
  653. static void partext()
  654. {
  655.     long i, l;
  656.     struct resbuf rbet, rbb;
  657.     struct resbuf *rb, *rent, *vb;
  658.     ads_name ename, vname;
  659.  
  660.     if (ptab != NULL) {
  661.         free(ptab);
  662.     }
  663.     ptab = NULL;
  664.     nptab = 0;
  665.  
  666.     /* Build the SSGET entity buffer chain to filter for block
  667.        insertions of the named block on the named layer. */
  668.  
  669.     rbet.restype = 0;                 /* Entity type */
  670.     rbet.resval.rstring = /*MSG0*/"INSERT";
  671.     rbet.rbnext = &rbb;
  672.     rbb.restype = 2;                  /* Block name */
  673.     rbb.resval.rstring = ParticleBlock;
  674.     rbb.rbnext = NULL;
  675.  
  676.     if (ads_ssget(/*MSG0*/"_X", NULL, NULL, &rbet, vname) != RTNORM) {
  677.         return;                       /* No definitions in database */
  678.     }
  679.  
  680.     /* We found one or more definitions.  Process  the  attributes
  681.        that  follow  each, plugging the values into material items
  682.        which get attached to the in-memory definition chain.  */
  683.  
  684.     if (ads_sslength(vname, &l) < 0)
  685.         l = 0;
  686.  
  687.     nptab = l;                        /* Save magnet count */
  688.     ntestp = (dualpend == 0.0) ? 1 : 2;  /* Number of test particles */
  689.     ptab = (particle *) alloc((nptab + ntestp) * sizeof(particle));
  690.     for (i = 0; i < l; i++) {
  691.         ads_name pname;
  692.  
  693.         ads_ssname(vname, i, ename);
  694.         memcpy(pname, ename, sizeof ename);
  695.  
  696.         memset(&pt, 0, sizeof pt);
  697.         while (True) {
  698.             ads_entnext(ename, ename);
  699.             rent = ads_entget(ename);
  700.             if (rent == NULL) {
  701.                 ads_printf(/*MSG30*/"PARTEXT: Can't read attribute.\n");
  702.                 break;
  703.             }
  704.             rb = entitem(rent, 0);    /* Entity type */
  705.             if (rb != NULL) {
  706.                 if (strcmp(rb->resval.rstring, /*MSG0*/"ATTRIB") != 0)
  707.                     break;
  708.                 rb = entitem(rent, 2);  /* Attribute tag */
  709.                 vb = entitem(rent, 1);  /* Attribute value */
  710.                 if (rb != NULL && vb != NULL) {
  711.                     if (strcmp(rb->resval.rstring, /*MSG31*/"PARTNAME") == 0) {
  712.                         strcpy(pt.partname, vb->resval.rstring);
  713.                     } else {
  714.                         int j;
  715.  
  716.                         for (j = 0; j < ELEMENTS(partatt); j++) {
  717.                             if (strcmp(rb->resval.rstring,
  718.                                        partatt[j].attname) == 0) {
  719.                                 *partatt[j].attvar = atof(vb->resval.rstring);
  720.                                 break;
  721.                             }
  722.                         }
  723.                     }
  724.                 }
  725.             }
  726.             ads_relrb(rent);
  727.         }
  728.         entinfo(pname, pt.partblock, pt.position, &pt.radius, &pt.colour);
  729.         memcpy(&(ptab[(int) i]), &pt, sizeof(particle));
  730.     }
  731.  
  732.     /* Initialise the test particle(s) to their initial state. */
  733.  
  734.     strcpy(ptab[nptab].partname, /*MSG32*/"Pendulum");
  735.     Spoint(ptab[nptab].position, 4.0, 4.0,
  736.            (PendLength + 1.0) - sqrt(PendLength * PendLength - 4.0 * 4.0));
  737.     Spoint(ptab[nptab].velocity, 0.0, 0.0, 0.0);
  738.     Spoint(ptab[nptab].acceleration, 0.0, 0.0, 0.0);
  739.     Cpoint(ptab[nptab].lastpos, ptab[nptab].position);
  740.     ptab[nptab].strength = 1.0;
  741.     ptab[nptab].radius = 0.25;
  742.     ptab[nptab].colour = 1;           /* Red */
  743.  
  744.     if (ntestp > 1) {
  745.         double mactemp = (4.0 * 4.0 + (4.0 + dualpend) * (4.0 + dualpend));
  746.  
  747.         strcpy(ptab[nptab + 1].partname, /*MSG33*/"Displaced pendulum");
  748.         Spoint(ptab[nptab + 1].position, 4.0, 4.0 + dualpend,
  749.                (PendLength + 1.0) - sqrt(PendLength * PendLength - mactemp));
  750.         Spoint(ptab[nptab + 1].velocity, 0.0, 0.0, 0.0);
  751.         Spoint(ptab[nptab + 1].acceleration, 0.0, 0.0, 0.0);
  752.         Cpoint(ptab[nptab + 1].lastpos, ptab[nptab].position);
  753.         ptab[nptab + 1].strength = 1.0;
  754.         ptab[nptab + 1].radius = 0.25;
  755.         ptab[nptab + 1].colour = 3;   /* Green */
  756.     }
  757.  
  758.     ads_ssfree(vname);                /* Release selection set */
  759. }
  760.  
  761. /*  SUMVEC  --  Add a linear multiple to another vector, a = b + t * c.  */
  762.  
  763. static void sumvec(ap, bp, t, cp)
  764.   ads_real *ap, *bp, *cp;
  765.   ads_real t;
  766. {
  767.     ap[X] = bp[X] + t * cp[X];
  768.     ap[Y] = bp[Y] + t * cp[Y];
  769.     ap[Z] = bp[Z] + t * cp[Z];
  770. }
  771.  
  772. /*  VARBLOCKDEF  --  Create the block that carries the attributes that
  773.                      define the modal variables. */
  774.  
  775. static void varblockdef()
  776. {
  777.     int i;
  778.     ads_name varbss;
  779.     ads_name ename;
  780.     ads_point atx;
  781.  
  782.     atx[X] = 4.0;
  783.     atx[Y] = 4.0;
  784.     atx[Z] = 0.0;
  785.     ads_ssadd(NULL, NULL, varbss);
  786.  
  787.     CommandB();
  788.     PushVar(/*MSG0*/"AFLAGS", saflags, 1, rint);  /* Invisible */
  789.  
  790.     for (i = 0; i < ELEMENTS(varatt); i++) {
  791.         atx[Y] -= 0.25;
  792.         ads_command(RTSTR, /*MSG0*/"_.ATTDEF",
  793.                     RTSTR, "", RTSTR, varatt[i].attname,
  794.                     RTSTR, varatt[i].attprompt, RTSTR, "0", RT3DPOINT, atx,
  795.                     RTREAL, 0.2, RTREAL, 0.0, RTNONE);
  796.         ads_entlast(ename);
  797.         ads_ssadd(ename, varbss, varbss);
  798.     }
  799.     PopVar(/*MSG0*/"AFLAGS", saflags);
  800.  
  801.     ads_command(RTSTR, /*MSG0*/"_.BLOCK", RTSTR, ModalBlock, RTSTR, "4,4",
  802.                 RTPICKS, varbss, RTSTR, "", RTNONE);
  803.     CommandE();
  804.     ads_ssfree(varbss);
  805. }
  806.  
  807. /*  SAVEMODES  --  Save application modes as attributes of the mode
  808.                    block.  */
  809.  
  810. static void savemodes()
  811. {
  812.     int i;
  813.     struct resbuf rbb;
  814.     ads_name ename, vname;
  815.     long l;
  816.  
  817.     /* Build the SSGET entity buffer chain to filter for block
  818.        insertions of the named block on the named layer. */
  819.  
  820.     rbb.restype = 2;                  /* Block name */
  821.     rbb.resval.rstring = ModalBlock;
  822.     rbb.rbnext = NULL;
  823.  
  824.     if (ads_ssget(/*MSG0*/"_X", NULL, NULL, &rbb, vname) == RTNORM) {
  825.  
  826.         /* Delete all definitions found. */
  827.  
  828.         if (ads_sslength(vname, &l) < 0)
  829.             l = 0;
  830.  
  831.         if (l > 0) {
  832.             CommandB();
  833.             ads_command(RTSTR, /*MSG0*/"_.ERASE",
  834.                         RTPICKS, vname, RTSTR, "", RTNONE);
  835.             CommandE();
  836.         }
  837.         ads_ssfree(vname);            /* Release selection set */
  838.     }
  839.  
  840.     /* Now insert the modal variable block, attaching the
  841.        mode variables to its attributes. */
  842.  
  843.     CommandB();
  844.     ads_command(RTSTR, /*MSG0*/"_.INSERT", RTSTR, ModalBlock,
  845.                 RTSTR, "0,0", RTSTR, "1", RTSTR, "1", RTSTR, "0",
  846.                 RTNONE);
  847.     for (i = 0; i < ELEMENTS(varatt); i++) {
  848.         char attval[20];
  849.  
  850. #define YorN(x) ((x) ? /*MSG34*/"Yes" : /*MSG35*/"No")
  851.         switch (varatt[i].attvari) {
  852.         case 1:
  853.             strcpy(attval, YorN(drawonly));
  854.             break;
  855.         case 2:
  856.             strcpy(attval, YorN(showtime));
  857.             break;
  858.         case 3:
  859.             strcpy(attval, YorN(showstep));
  860.             break;
  861.         case 4:
  862.             sprintf(attval, "%.12g", stepsize);
  863.             break;
  864.         case 5:
  865.             sprintf(attval, "%.12g", dualpend);
  866.             break;
  867.         }
  868. #undef YorN
  869.         ads_command(RTSTR, attval, RTNONE);
  870.     }
  871.     ads_entlast(ename);
  872.     ads_command(RTSTR, /*MSG0*/"_.CHPROP", RTENAME, ename, RTSTR, "",
  873.                 RTSTR, /*MSG0*/"_layer", RTSTR, FrozenLayer, RTSTR, "",
  874.                 RTNONE);
  875.     CommandE();
  876. }
  877.  
  878. /*  BOOLVAR  --  Determine Boolean value from attribute string.  We
  879.                  recognise numbers (nonzero means True) and the
  880.                  strings "True", "False", "Yes", and "No",
  881.                  in either upper or lower case.  */
  882.  
  883. static Boolean boolvar(varname, s)
  884.   char *varname, *s;
  885. {
  886.     char ch = *s;
  887.  
  888.     if (islower(ch))
  889.         ch = toupper(ch);
  890.     if (isdigit(ch)) {
  891.         return (atoi(s) != 0) ? True : False;
  892.     }
  893.  
  894.     switch (ch) {
  895.     case /*MSG36*/'Y':
  896.     case /*MSG37*/'T':
  897.         return True;
  898.     case /*MSG38*/'N':
  899.     case /*MSG39*/'F':
  900.         return False;
  901.     default:
  902.         ads_printf(/*MSG40*/"Invalid Boolean value for %s: %s.\n", varname, s);
  903.         break;
  904.     }
  905.     return False;
  906. }
  907.  
  908. /*  VARSET  --  Set modal variables from the block attributes.  */
  909.  
  910. static void varset()
  911. {
  912.     struct resbuf rbb;
  913.     ads_name vname;
  914.     long l;
  915.  
  916.     /* Build the SSGET entity buffer chain to filter for block
  917.        insertions of the named block on the named layer. */
  918.  
  919.     rbb.restype = 2;                  /* Block name */
  920.     rbb.resval.rstring = ModalBlock;
  921.     rbb.rbnext = NULL;
  922.  
  923.     if (ads_ssget(/*MSG0*/"_X", NULL, NULL, &rbb, vname) == RTNORM) {
  924.  
  925.         if (ads_sslength(vname, &l) < 0)
  926.             l = 0;
  927.  
  928.         if (l > 0) {
  929.             ads_name ename;
  930.  
  931.             if (ads_ssname(vname, 0L, ename) == RTNORM) {
  932.                 while (ads_entnext(ename, ename) == RTNORM) {
  933.                     int i;
  934.                     struct resbuf *rb = ads_entget(ename),
  935.                                   *et, *at, *av;
  936.  
  937.                     if (rb == NULL)
  938.                         break;
  939.                     et = entitem(rb, 0);
  940.                     assert(et != NULL);
  941.                     if (strcmp(et->resval.rstring, /*MSG0*/"ATTRIB") == 0) {
  942.                         at = entitem(rb, 2);  /* Attribute tag */
  943.                         av = entitem(rb, 1);  /* Attribute value */
  944.                         if (at == NULL || av == NULL)
  945.                             break;
  946.                         for (i = 0; i < ELEMENTS(varatt); i++) {
  947.                             if (strcmp(at->resval.rstring,
  948.                                        varatt[i].attname) == 0) {
  949.                                 switch (varatt[i].attvari) {
  950.                                 case 1:
  951.                                     drawonly = boolvar(at->resval.rstring,
  952.                                                        av->resval.rstring);
  953.                                     break;
  954.                                 case 2:
  955.                                     showtime = boolvar(at->resval.rstring,
  956.                                                        av->resval.rstring);
  957.                                     break;
  958.                                 case 3:
  959.                                     showstep = boolvar(at->resval.rstring,
  960.                                                        av->resval.rstring);
  961.                                     break;
  962.                                 case 4:
  963.                                     stepsize = 0.1;    /* Default if error */
  964.                                     sscanf(av->resval.rstring, "%lf",
  965.                                            &stepsize);
  966.                                     break;
  967.                                 case 5:
  968.                                     dualpend = 0.0;    /* Default if error */
  969.                                     sscanf(av->resval.rstring, "%lf",
  970.                                            &dualpend);
  971.                                     break;
  972.                                 }
  973.                             }
  974.                         }
  975.                     } else if (strcmp(et->resval.rstring,
  976.                                       /*MSG0*/"SEQEND") == 0) {
  977.                         ads_relrb(rb);
  978.                         break;
  979.                     }
  980.                     ads_relrb(rb);
  981.                 }
  982.             }
  983.         }
  984.         ads_ssfree(vname);            /* Release selection set */
  985.     }
  986. }
  987.  
  988. /*  INITACAD  --  Initialise the required modes in the AutoCAD
  989.                   drawing.  */
  990.  
  991. static Boolean initacad(reset)
  992.   Boolean reset;
  993. {
  994.     static Boolean initdone, initok;
  995.  
  996.     if (reset) {
  997.         initdone = False;
  998.         initok = True;
  999.     } else {
  1000.         if (!initdone) {
  1001.             struct resbuf rb;
  1002.             struct resbuf *ep;
  1003.             ads_point pzero;
  1004.  
  1005.             initok = False;
  1006.  
  1007.             /* Reset the program modes to standard values upon
  1008.                entry to the drawing editor. */
  1009.  
  1010.             stepno = 0;               /* Step 0 */
  1011.             numsteps = 1000;          /* Default number of steps */
  1012.             simtime = 0.0;            /* No elapsed time */
  1013.             stepsize = 0.1;           /* Default step size */
  1014.             dualpend = 0.0;           /* Default 2nd pendulum offset: none */
  1015.             drawonly = True;          /* Draw using ads_grdraw() */
  1016.             showstep = True;          /* Show step number */
  1017.             showtime = True;          /* Show elapsed time */
  1018.  
  1019.             /* Enable  handles  if  they aren't  already on.  We use
  1020.                handles to link  from  our  in-memory  table  to  the
  1021.                magnets in  the  AutoCAD  database,  and we also keep
  1022.                track of the reference frame entity  by  its  handle. */
  1023.  
  1024.             ads_getvar(/*MSG0*/"HANDLES", &rb);
  1025.             if (!rb.resval.rint) {
  1026.                 CommandB();
  1027.                 ads_command(RTSTR, /*MSG0*/"_.handles",
  1028.                             RTSTR, /*MSG0*/"_on", RTNONE);
  1029.                 CommandE();
  1030.                 ads_getvar(/*MSG0*/"HANDLES", &rb);
  1031.                 if (!rb.resval.rint) {
  1032.                     ads_printf(/*MSG41*/"Cannot enable handles.\n");
  1033.                     initdone = True;
  1034.                     return (initok = False);
  1035.                 }
  1036.             }
  1037.  
  1038.             /* Create required drawing objects. */
  1039.  
  1040.             /* If there isn't already a "frozen solid" layer, create one. */
  1041.  
  1042.             if ((ep = ads_tblsearch(/*MSG0*/"LAYER",
  1043.                                     FrozenLayer, False)) == NULL) {
  1044.                 CommandB();
  1045.                 ads_command(RTSTR, /*MSG0*/"_.LAYER",
  1046.                             RTSTR, /*MSG0*/"_NEW", RTSTR, FrozenLayer,
  1047.                             RTSTR, /*MSG0*/"_FREEZE",
  1048.                             RTSTR, FrozenLayer, RTSTR, "", RTNONE);
  1049.                 ads_command(RTSTR, /*MSG0*/"_.LAYER",
  1050.                             RTSTR, /*MSG0*/"_NEW", RTSTR, OrbitLayer,
  1051.                             RTSTR, "", RTNONE);
  1052.                 CommandE();
  1053.                 defmagblk();          /* Define the particle block */
  1054.             } else {
  1055.                 ads_relrb(ep);
  1056.             }
  1057.  
  1058.             /* Create the block definition for our modal variables. */
  1059.  
  1060.             if ((ep = ads_tblsearch(/*MSG0*/"BLOCK",
  1061.                                     ModalBlock, False)) == NULL) {
  1062.                 varblockdef();
  1063.                 savemodes();
  1064.             } else {
  1065.                 ads_relrb(ep);
  1066.                 varset();             /* Load modals from block */
  1067.             }
  1068.  
  1069.             /* Establish default initial view. */
  1070.  
  1071.             CommandB();
  1072.             Spoint(pzero, 0.0, 0.0, 0.0);
  1073.             ads_command(RTSTR, /*MSG0*/"_.ZOOM",
  1074.                         RTSTR, /*MSG0*/"_C", RT3DPOINT, pzero,
  1075.                         RTREAL, 2.0 * (PendLength / PcX), RTNONE);
  1076.             CommandE();
  1077.  
  1078.             initdone = initok = True;
  1079.         }
  1080.     }
  1081.     return initok;
  1082. }
  1083.  
  1084. /*  ADDMAG  --  Insert magnet block with attributes.  */
  1085.  
  1086. static void addmag(mname, pos, strength)
  1087.   char *mname;
  1088.   ads_point pos;
  1089.   ads_real strength;
  1090. {
  1091.     ads_real size;
  1092.     char smas[32];
  1093.  
  1094.     size = 1.0;
  1095.     CommandB();
  1096.     sprintf(smas, "%.12g", strength);
  1097.     ads_command(RTSTR, /*MSG0*/"_.INSERT", RTSTR, ParticleBlock, RTPOINT, pos,
  1098.                 RTREAL, size, RTREAL, size, RTREAL, 0.0,
  1099.                 RTSTR, mname,
  1100.                 RTSTR, smas, RTNONE);
  1101.     CommandE();
  1102. }
  1103.  
  1104. /*  MAGNET  --  Add magnet to database.  */
  1105.  
  1106. static void magnet()
  1107. {
  1108.     static int magseq = 0;
  1109.     ads_point pos;
  1110.     ads_real strength;
  1111.     char pname[134];
  1112.  
  1113.     sprintf(pname, /*MSG42*/"Magnet %d", ++magseq);
  1114.     ads_initget(1 + 8 + 16, NULL);
  1115.     if (ads_getpoint(NULL, /*MSG43*/"\nPosition: ", pos) != RTNORM)
  1116.         return;
  1117.     ads_initget(2, NULL);             /* No zero-strength magnets */
  1118.     switch (ads_getreal(/*MSG44*/"\nStrength of magnet <1>: ", &strength)) {
  1119.     case RTNONE:
  1120.         strength = 1.0;
  1121.     case RTNORM:
  1122.         break;
  1123.     default:
  1124.         return;
  1125.     }
  1126.     addmag(pname, pos, strength);
  1127.     ready++;
  1128. }
  1129.  
  1130. /*  DEMO  --  Load canned demo case.  */
  1131.  
  1132. static void demo()
  1133. {
  1134.     int i, j;
  1135.  
  1136.     typedef struct {                  /* Demo magnet table item */
  1137.        char *dmname;                  /* Magnet name */
  1138.        ads_point dmpos;               /* Position */
  1139.        ads_real dmstrength;           /* Magnetic strength */
  1140.        int dmcolour;                  /* Colour */
  1141.     } dmag;
  1142.  
  1143.     typedef struct {
  1144.        char *dname;                   /* Name of demo case */
  1145.        ads_real dnstep;               /* Default number of steps */
  1146.        ads_real dssize;               /* Step size */
  1147.        ads_real dualp;                /* Dual pendulum offset */
  1148.        dmag *dmtab;                   /* Table of magnets */
  1149.        int dmtabl;                    /* Number of magnets in table */
  1150.     } demodesc;
  1151.  
  1152.     static dmag interlom[] = {
  1153.        {/*MSG45*/"Demo magnet 1", {3, 3, 0}, 1, 2},
  1154.        {/*MSG46*/"Demo magnet 2", {0, 3, 0}, 2, 3},
  1155.        {/*MSG47*/"Demo magnet 3", {-1, -2, 0}, 1, 1},
  1156.        {/*MSG48*/"Demo magnet 3", {2, -3, 0}, 1, 5}
  1157.       };
  1158.     static demodesc dmt[] = {
  1159.        {/*MSG49*/"Capt. Magneto",
  1160.         3000, 0.1, 0.001,
  1161.         interlom, ELEMENTS(interlom)}
  1162.       };
  1163.  
  1164.     partext();
  1165.     if (nptab > 0) {
  1166.         ads_printf(/*MSG50*/"\n\
  1167. Demo can be created only in an empty drawing.\n");
  1168.         return;
  1169.     }
  1170.  
  1171. #ifndef MULTIPLE_DEMOS
  1172.     i = 0;                            /* Always use first demo */
  1173. #else
  1174.     ads_textscr();                    /* Flip to text screen */
  1175.     ads_printf(/*MSG51*/"Available demonstration cases:\n");
  1176.     for (i = 0; i < ELEMENTS(dmt); i++) {
  1177.         ads_printf("\n%d.  %s", i + 1, dmt[i].dname);
  1178.     }
  1179.     ads_printf("\n");
  1180.     ads_initget(2 + 4, NULL);
  1181.     switch (ads_getint(/*MSG52*/"\nWhich demonstration? ", &i)) {
  1182.     case RTNORM:
  1183.         i--;
  1184.         if (i < ELEMENTS(dmt))
  1185.             break;
  1186.     default:
  1187.         return;
  1188.     }
  1189. #endif
  1190.  
  1191.     for (j = 0; j < dmt[i].dmtabl; j++) {
  1192.         CommandB();
  1193.         ads_command(RTSTR, /*MSG0*/"_.COLOUR",
  1194.                     RTSHORT, dmt[i].dmtab[j].dmcolour, RTNONE);
  1195.         CommandE();
  1196.         addmag(dmt[i].dmtab[j].dmname, dmt[i].dmtab[j].dmpos,
  1197.                dmt[i].dmtab[j].dmstrength);
  1198.     }
  1199.  
  1200.     numsteps = dmt[i].dnstep;
  1201.     stepsize = dmt[i].dssize;
  1202.     dualpend = dmt[i].dualp;
  1203.     savemodes();
  1204.     CommandB();
  1205.     ads_command(RTSTR, /*MSG0*/"_.COLOUR", RTSTR, /*MSG0*/"_BYLAYER", RTNONE);
  1206.     CommandE();
  1207.     partext();
  1208.     ready++;
  1209. }
  1210.  
  1211. /*  MAGEDIT  --  Edit properties of an existing magnet.  */
  1212.  
  1213. static void magedit()
  1214. {
  1215.     int i;
  1216.     ads_name ename;
  1217.     ads_point ppt;
  1218.  
  1219.     partext();                        /* Update in-memory magnet table */
  1220.     i = nptab; 
  1221.     if (ads_entsel(/*MSG53*/"Select magnet to edit: ", ename, ppt) == RTNORM) {
  1222.         struct resbuf *eb = ads_entget(ename), *ha;
  1223.  
  1224.         if ((ha = entitem(eb, 5)) != NULL) {
  1225.             for (i = 0; i < nptab; i++) {
  1226.                 if (strcmp(ptab[i].partblock, ha->resval.rstring) == 0) {
  1227.                     break;
  1228.                 }
  1229.             }
  1230.         }
  1231.     }
  1232.     if (i < nptab) {
  1233.         if (stepno > 0) {
  1234.             reset();
  1235.         }
  1236.         CommandB();
  1237.         ads_command(RTSTR, /*MSG0*/"_.DDATTE", RTENAME, ename, RTNONE);
  1238.         CommandE();
  1239.         partext();
  1240.     } else {
  1241.         ads_printf(/*MSG54*/"\nNo magnet selected.\n");
  1242.     }
  1243. }
  1244.  
  1245. /*  RESET  --  Reset simulation, erasing all pendulum paths.  */
  1246.  
  1247. void reset()
  1248. {
  1249.     struct resbuf rbb;
  1250.     ads_name vname;
  1251.     long l;
  1252.  
  1253.     /* Build the SSGET entity buffer chain to select all
  1254.        entities on the orbital path layer. */
  1255.  
  1256.     rbb.restype = 8;                  /* Layer name */
  1257.     rbb.resval.rstring = OrbitLayer;
  1258.     rbb.rbnext = NULL;
  1259.  
  1260.     if (ads_ssget(/*MSG0*/"_X", NULL, NULL, &rbb, vname) == RTNORM) {
  1261.  
  1262.         /* We found one or more definitions.  Delete them. */
  1263.  
  1264.         if (ads_sslength(vname, &l) < 0)
  1265.             l = 0;
  1266.  
  1267.         if (l > 0) {
  1268.             CommandB();
  1269.             ads_command(RTSTR, /*MSG0*/"_.ERASE",
  1270.                         RTPICKS, vname, RTSTR, "", RTNONE);
  1271.             CommandE();
  1272.         }
  1273.         ads_ssfree(vname);            /* Release selection set */
  1274.     }
  1275.     if (drawonly) {
  1276.         ads_grtext(-3, NULL, 0);      /* Clear time display */
  1277.         ads_redraw(NULL, 0);
  1278.     }
  1279.     stepno = 0;                       /* Reset step number */
  1280. }
  1281.  
  1282. /*  RUN  --  Run simulation.  This can be called as an interactive
  1283.              command, RUN, or functionally as (C:RUN <steps/-time>).  */
  1284.  
  1285. static void run()
  1286. {
  1287.     int i, j, k, n = 0, lastcolour = -1;
  1288.     Boolean resume = (stepno > 0) ? True : False, timelimit;
  1289.     ads_real r, endtime = 0;
  1290.     char rpt[80];
  1291.     struct resbuf clayer, ocolour;
  1292.     struct resbuf *rb = ads_getargs();
  1293.  
  1294.     if (!ready) {
  1295.         ads_printf( "\nSet up some magnets first or run Demo.\n");
  1296.         return;
  1297.     }
  1298.         
  1299.     /* If an argument was supplied, set the simulation length
  1300.        from it rather than asking the user interactively. */
  1301.  
  1302.     if (rb != NULL) {
  1303.         Boolean argok = True;
  1304.  
  1305.         switch (rb->restype) {
  1306.         case RTREAL:
  1307.             numsteps = rb->resval.rreal;
  1308.             break;
  1309.         case RTSHORT:
  1310.             numsteps = rb->resval.rint;
  1311.             break;
  1312.         default:
  1313.             ads_fail(/*MSG55*/"RUN: Improper argument type");
  1314.             argok = False;
  1315.             break;
  1316.         }
  1317.         if (argok) {
  1318.             if (numsteps == 0) {
  1319.                 ads_fail(/*MSG56*/"RUN: Argument must be nonzero");
  1320.                 argok = False;
  1321.             } else {
  1322.                 if (rb->rbnext != NULL) {
  1323.                     ads_fail(/*MSG57*/"RUN: Too many arguments");
  1324.                     argok = False;
  1325.                 }
  1326.             }
  1327.         }
  1328.         if (!argok)
  1329.             return;
  1330.         functional = True;            /* Mark called as a function */
  1331.     } else {
  1332.  
  1333.         /* Ask user for length of simulation in years or number of
  1334.            integration steps to perform. */
  1335.  
  1336.         sprintf(rpt, /*MSG58*/"Steps or (-simulated seconds) <%.12g>: ",
  1337.                 numsteps);
  1338.         ads_initget(2, NULL);
  1339.         switch (ads_getreal(rpt, &r)) {
  1340.         case RTCAN:                   /* Control C */
  1341.             return;
  1342.         case RTNONE:                  /* Null input */
  1343.             break;
  1344.         case RTNORM:                  /* Number entered */
  1345.             numsteps = r;
  1346.             break;
  1347.         }
  1348.     }
  1349.  
  1350.     /* If we're starting a new simulation rather than resuming one
  1351.        already underway, reset  the  simulated  time  counter  and
  1352.        extract the current particle information from the database. */
  1353.  
  1354.     if (!resume) {
  1355.         ads_point bpx;
  1356.  
  1357.         Spoint(bpx, 0.0, 0.0, 0.0);
  1358.         simtime = 0.0;
  1359.         partext();
  1360.         ads_initget(1 + 8, NULL);
  1361.         if (ads_getpoint(bpx, /*MSG59*/"Pendulum start position: ", 
  1362.                 ptab[nptab].position) < 0)
  1363.         {
  1364.             Spoint(ptab[nptab].position, 4.0, 4.0, (PendLength + 1.0) -
  1365.                                                    sqrt(PendLength *
  1366.                                                         PendLength -
  1367.                                                         4.0 * 4.0));
  1368.         } else {
  1369.             ptab[nptab].position[Z] = (PendLength + 1.0) -
  1370.                                       sqrt(PendLength * PendLength -
  1371.                                            (ptab[nptab].position[X] *
  1372.                                             ptab[nptab].position[X] +
  1373.                                             ptab[nptab].position[Y] *
  1374.                                             ptab[nptab].position[Y]));
  1375.         }
  1376.         Cpoint(ptab[nptab].lastpos, ptab[nptab].position);
  1377.         if (ntestp > 1) {
  1378.             Cpoint(ptab[nptab + 1].position, ptab[nptab].position);
  1379.             ptab[nptab + 1].position[Y] += dualpend;
  1380.             ptab[nptab + 1].position[Z] = (PendLength + 1.0) -
  1381.                                           sqrt(PendLength * PendLength -
  1382.                                                (ptab[nptab + 1].position[X] *
  1383.                                                 ptab[nptab].position[X] +
  1384.                                                 ptab[nptab + 1].position[Y] *
  1385.                                                 ptab[nptab].position[Y]));
  1386.             Cpoint(ptab[nptab + 1].lastpos, ptab[nptab + 1].position);
  1387.         }
  1388.     }
  1389.  
  1390.     if (numsteps > 0) {
  1391.         timelimit = False;
  1392.         n = numsteps;
  1393.     } else {
  1394.         timelimit = True;
  1395.         endtime = simtime - numsteps;
  1396.     }
  1397.  
  1398.     /* Save original layer and colour */
  1399.  
  1400.     ads_getvar(/*MSG0*/"CLAYER", &clayer);
  1401.     ads_getvar(/*MSG0*/"CECOLOR", &ocolour);
  1402.     /* AutoCAD  currently returns  "human readable" colour strings
  1403.        like "1 (red)" for the standard colours.  Trim  the  string
  1404.        at  the  first space to guarantee we have a valid string to
  1405.        restore the colour later.  */
  1406.     if (strchr(ocolour.resval.rstring, ' ') != NULL)
  1407.         *strchr(ocolour.resval.rstring, ' ') = EOS;
  1408.  
  1409.     CommandB();
  1410.     ads_command(RTSTR, /*MSG0*/"_.UCS", RTSTR, /*MSG0*/"_World", RTNONE);
  1411.     if (!drawonly) {
  1412.         ads_command(RTSTR, /*MSG0*/"_.LAYER",
  1413.                     RTSTR, /*MSG0*/"_SET", RTSTR, OrbitLayer,
  1414.                     RTSTR, "", RTNONE);
  1415.     }
  1416.  
  1417.     /* Initialise last plotted position for each particle. */
  1418.  
  1419.     for (i = 0; i < nptab; i++) {
  1420.         Cpoint(ptab[i].lastpos, ptab[i].position);
  1421.     }
  1422.  
  1423.     while ((timelimit ? (simtime < endtime) : (n-- > 0)) &&
  1424.            !ads_usrbrk()) {
  1425.  
  1426.         /* Update acceleration for all test bodies. */
  1427.  
  1428.         for (i = nptab; i < (nptab + ntestp); i++) {
  1429.             ads_real dist, pacce;
  1430.             ads_point pzero;
  1431.  
  1432.             /* Initialise acceleration to that resulting from the
  1433.                restoring force exerted by the current position of
  1434.                the pendulum. */
  1435.  
  1436.             Spoint(pzero, 0.0, 0.0, 0.0);
  1437.             dist = ads_distance(ptab[i].position, pzero);
  1438.             pacce = - ptab[i].strength * (dist / PendLength);
  1439.             Spoint(ptab[i].acceleration,
  1440.                    pacce * (ptab[i].position[X] / dist),
  1441.                    pacce * (ptab[i].position[Y] / dist),
  1442.                    0.0);
  1443.             for (j = 0; j < nptab; j++) {
  1444.                 if (i != j) {
  1445.                     ads_real vdist = ads_distance(ptab[i].position,
  1446.                                                   ptab[j].position),
  1447.                              /* F = K (m1 m2) / r^2 */
  1448.                              force = ((ptab[i].strength * ptab[j].strength) /
  1449.                                       (vdist * vdist)),
  1450.                              /* F = ma, hence a = F/m.  We take the mass to
  1451.                                 be 1, so the acceleration is just equal to
  1452.                                 the force. */
  1453.                              accel = force;
  1454.  
  1455.                     if (vdist == 0.0) {
  1456.                         ads_printf(/*MSG61*/"Collision between %s and %s\n",
  1457.                                    ptab[i].partname, ptab[j].partname);
  1458.                     } else {
  1459.                         /* Update vector components of acceleration. */
  1460.                         for (k = X; k <= Z; k++) {
  1461.                             ptab[i].acceleration[k] +=
  1462.                                accel * (ptab[i].position[k] -
  1463.                                         ptab[j].position[k]) / vdist;
  1464.                         }
  1465.                     }
  1466.                 }
  1467.             }
  1468.         }
  1469.  
  1470.         /* Update velocity for all bodies. */
  1471.  
  1472.         for (i = nptab; i < (nptab + ntestp); i++) {
  1473.             sumvec(ptab[i].velocity, ptab[i].velocity,
  1474.                    stepsize, ptab[i].acceleration);
  1475.         }
  1476.  
  1477.         /* Update position for all bodies. */
  1478.  
  1479.         for (i = nptab; i < (nptab + ntestp); i++) {
  1480.             ads_real pdx;
  1481.  
  1482.             sumvec(ptab[i].position, ptab[i].position,
  1483.                    stepsize, ptab[i].velocity);
  1484.             pdx = ptab[i].position[X] * ptab[i].position[X] +
  1485.                   ptab[i].position[Y] * ptab[i].position[Y];
  1486.             if (pdx > (PendLength * PendLength)) {
  1487.                 ads_real patx;
  1488.  
  1489.                 ads_printf(/*MSG62*/"\
  1490. Boink!!  Hit the stop--magnets too strong!\n");
  1491.                 Spoint(ptab[i].acceleration, 0.0, 0.0, 0.0);
  1492.                 Spoint(ptab[i].velocity, 0.0, 0.0, 0.0);
  1493.                 patx = atan2(ptab[i].position[Y], ptab[i].position[X]);
  1494.                 ptab[i].position[X] = PendLength * cos(patx);
  1495.                 ptab[i].position[Y] = PendLength * sin(patx);
  1496.                 pdx = PendLength * PendLength;
  1497.             }
  1498.             ptab[i].position[Z] = (PendLength + 1.0) -
  1499.                                   sqrt(PendLength * PendLength - pdx);
  1500.         }
  1501.  
  1502.         /* Display motion since last update. */
  1503.  
  1504.         for (i = nptab; i < (nptab + ntestp); i++) {
  1505.             if (drawonly) {
  1506.                 ads_grdraw(ptab[i].lastpos, ptab[i].position,
  1507.                            ptab[i].colour, False);
  1508.             } else {
  1509.                 if (lastcolour != ptab[i].colour) {
  1510.                     ads_command(RTSTR, /*MSG0*/"_.COLOUR",
  1511.                                 RTSHORT, ptab[i].colour,
  1512.                                 RTNONE);
  1513.                     lastcolour = ptab[i].colour;
  1514.                 }
  1515.                 ads_command(RTSTR, /*MSG0*/"_.LINE", RTPOINT, ptab[i].lastpos,
  1516.                             RTPOINT, ptab[i].position, RTSTR, "", RTNONE);
  1517.             }
  1518.             Cpoint(ptab[i].lastpos, ptab[i].position);
  1519.         }
  1520.  
  1521.         /* Update the step number, simulated time, and display
  1522.            the values on the status line if selected. */
  1523.  
  1524.         stepno++;
  1525.         simtime += stepsize;
  1526.         rpt[0] = EOS;
  1527.         if (showtime) {
  1528.             sprintf(rpt, /*MSG63*/"Time %.2f", simtime);
  1529.         }
  1530.         if (showstep) {
  1531.             if (rpt[0] != EOS)
  1532.                 strcat(rpt, "  ");
  1533.             sprintf(rpt + strlen(rpt), /*MSG64*/"Step %ld", stepno);
  1534.         }
  1535.         if (rpt[0] != EOS)
  1536.             ads_grtext(-2, rpt, False);
  1537.     }
  1538.  
  1539.     /* Restore original layer, colour, and UCS. */
  1540.  
  1541.     if (!drawonly) {
  1542.         ads_command(RTSTR, /*MSG0*/"_.LAYER", RTSTR, /*MSG0*/"_SET",
  1543.                     RTSTR, clayer.resval.rstring, RTSTR, "", RTNONE);
  1544.         ads_command(RTSTR, /*MSG0*/"_.COLOUR",
  1545.                     RTSTR, ocolour.resval.rstring, RTNONE);
  1546.     }
  1547.     free(clayer.resval.rstring);
  1548.     free(ocolour.resval.rstring);
  1549.     ads_command(RTSTR, /*MSG0*/"_.UCS", RTSTR, /*MSG0*/"_Prev", RTNONE);
  1550.     CommandE();
  1551.  
  1552.     /* If called as a function, return the simulation end time.
  1553.        If called as a command, print the end time and step count. */
  1554.  
  1555.     if (functional)
  1556.         ads_retreal(simtime);
  1557.     else
  1558.         ads_printf(/*MSG65*/"\nEnd at step %ld, %.2f seconds.\n",
  1559.                    stepno, simtime);
  1560. }
  1561.  
  1562. /*  SETMAG  --  Set modal variables.  */
  1563.  
  1564. static void setmag()
  1565. {
  1566.     struct resbuf rbb;
  1567.     ads_name vname;
  1568.     long l;
  1569.  
  1570.     /* Build the SSGET entity buffer chain to filter for block
  1571.        insertions of the named block on the named layer. */
  1572.  
  1573.     rbb.restype = 2;                  /* Block name */
  1574.     rbb.resval.rstring = ModalBlock;
  1575.     rbb.rbnext = NULL;
  1576.  
  1577.     if (ads_ssget(/*MSG0*/"_X", NULL, NULL, &rbb, vname) == RTNORM) {
  1578.  
  1579.         /* Delete all definitions found. */
  1580.  
  1581.         if (ads_sslength(vname, &l) < 0)
  1582.             l = 0;
  1583.  
  1584.         if (l > 0) {
  1585.             ads_name ename;
  1586.  
  1587.             if (ads_ssname(vname, 0L, ename) == RTNORM) {
  1588.                 CommandB();
  1589.                 ads_command(RTSTR, /*MSG0*/"__DDATTE", RTENAME, ename, RTNONE);
  1590.                 CommandE();
  1591.                 varset();
  1592.             }
  1593.         }
  1594.         ads_ssfree(vname);            /* Release selection set */
  1595.     }
  1596. }
  1597.  
  1598.  
  1599. #ifdef  HIGHC
  1600.  
  1601. /*  Earlier versions of High C put abort() in the same module with exit();
  1602.     ADS defines its own exit(), so we have to define our own abort(): */
  1603.  
  1604. static void abort(void)
  1605. {
  1606.     ads_abort("");
  1607. }
  1608.  
  1609. #endif  /* HIGHC */
  1610.